macula densa subtype analysis
jeff karnsomprot
2025-04-18
1 Macula Densa Project
1.1 Objectives
Three Main Goals of this File
Produce Cleaner looking code.
Identify the amount of clusters there are
Identify the top genes
expressed in each of the clusters
2 Loading in Data sets + Library packages.
options(future.globals.maxSize = 74 * 1024^3) # 55 GB
getOption("future.globals.maxSize") #59055800320## [1] 794568949763 Analyzing DATASET
SO4@meta.data$sample <- factor(
SO4@meta.data$sample,
levels = c("SO1","SO4","SO2","SO3"))
DimPlot(SO4,split.by ="sample",group.by = "sample")
3.1 Viewing for different DEGs
markers.to.plot1 <- c(
"Atf3", #
"Egr1", #
"Fos", #
"Jun", #
"Junb", #
"Pappa2", #
"Cxcl10", #
"Cldn19", #
"Krt7", #
"Egf", #
"Ptger3",
"Ckb",
"Mcub",
"Fabp3",
"Foxq1",
"Vash2",
"Pamr1",
"Vegfa"
)
DotPlot(SO4,
features = markers.to.plot1,
dot.scale = 8,
dot.min = 0,
scale = FALSE,
scale.max = 100,
scale.min = 0,
col.min = -2.5,
group.by = "seurat_clusters",
col.max = 2.5)+
coord_flip()
3.2 Finding Markers for each MD Subtype + Saving Excel File
Idents(SO4) <- "supergroup"
MD_DEGs <- FindAllMarkers(SO4,
only.pos = TRUE,
logfc.threshold = 0.1,
min.pct = 0.1,
min.diff.pct = 0.1,
return.thresh = 0.05)## Calculating cluster type_1## Warning: The `slot` argument of `GetAssayData()` is deprecated as of SeuratObject 5.0.0.
## ℹ Please use the `layer` argument instead.
## ℹ The deprecated feature was likely used in the Seurat package.
## Please report the issue at <https://github.com/satijalab/seurat/issues>.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.## Warning: `PackageCheck()` was deprecated in SeuratObject 5.0.0.
## ℹ Please use `rlang::check_installed()` instead.
## ℹ The deprecated feature was likely used in the Seurat package.
## Please report the issue at <https://github.com/satijalab/seurat/issues>.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.## For a (much!) faster implementation of the Wilcoxon Rank Sum Test,
## (default method for FindMarkers) please install the presto package
## --------------------------------------------
## install.packages('devtools')
## devtools::install_github('immunogenomics/presto')
## --------------------------------------------
## After installation of presto, Seurat will automatically use the more
## efficient implementation (no further action necessary).
## This message will be shown once per session## Calculating cluster type_2## Calculating cluster type_3## Calculating cluster type_4## Calculating cluster type_5## Calculating cluster type_6MD_DEGs <- MD_DEGs %>%
arrange(desc(avg_log2FC)) %>%
dplyr::select(gene, everything())
# Split by cluster (ident column)
marker_list <- split(MD_DEGs, MD_DEGs$cluster)
# Sort names alphanumerically
sorted_cluster_names <- names(marker_list)[order(as.numeric(names(marker_list)))]## Warning in eval(quote(list(...)), env): NAs introduced by coercion# Create a workbook
wb <- createWorkbook()
# Add each cluster as a new worksheet
for (cluster_name in sorted_cluster_names) {
addWorksheet(wb, sheetName = cluster_name)
writeData(wb, sheet = cluster_name, marker_list[[cluster_name]])
}
date <- format(Sys.Date(), "%Y%m%d")
# Save workbook
saveWorkbook(wb, here("jk_code", paste0(date, "_", "_SO4_FindAllMarkers_By_SubType2.xlsx")), overwrite = TRUE)
MD_DEGs %>%
group_by(cluster) %>%
dplyr::filter(avg_log2FC > 1, p_val_adj < 0.05) %>%
slice_head(n = 10) %>%
ungroup() -> top10
# save as excel file
DoHeatmap(SO4, features = top10$gene) + NoLegend()## Warning in DoHeatmap(SO4, features = top10$gene): The following features were
## omitted as they were not found in the scale.data slot for the SCT assay: H2-Q6,
## Gbp3, Ifi47, Tmem265
3.3 Finding Markers for each MD Subtype + Saving Excel File
Idents(SO4) <- "manual_groups"
MD_DEGs2 <- FindAllMarkers(SO4,
only.pos = TRUE,
logfc.threshold = 0.1,
min.pct = 0.1,
min.diff.pct = 0.1,
return.thresh = 0.05)## Calculating cluster type_1a## Calculating cluster type_1b## Calculating cluster type_1c## Calculating cluster type_2a## Calculating cluster type_2b## Calculating cluster type_3a## Calculating cluster type_3b## Calculating cluster type_4## Calculating cluster type_5## Calculating cluster type_6MD_DEGs2 <- MD_DEGs2 %>%
arrange(desc(avg_log2FC)) %>%
dplyr::select(gene, everything())
# Split by cluster (ident column)
marker_list <- split(MD_DEGs2, MD_DEGs2$cluster)
# Sort names alphanumerically
sorted_cluster_names <- names(marker_list)[order(as.numeric(names(marker_list)))]## Warning in eval(quote(list(...)), env): NAs introduced by coercion# Create a workbook
wb <- createWorkbook()
# Add each cluster as a new worksheet
for (cluster_name in sorted_cluster_names) {
addWorksheet(wb, sheetName = cluster_name)
writeData(wb, sheet = cluster_name, marker_list[[cluster_name]])
}
date <- format(Sys.Date(), "%Y%m%d")
# Save workbook
saveWorkbook(wb, here("jk_code", paste0(date, "_", "_SO4_FindAllMarkers_By_Sub_SubType2.xlsx")), overwrite = TRUE)
MD_DEGs %>%
group_by(cluster) %>%
dplyr::filter(avg_log2FC > 1, p_val_adj < 0.05) %>%
slice_head(n = 10) %>%
ungroup() -> top10
# save as excel file
DoHeatmap(SO4, features = top10$gene) + NoLegend()## Warning in DoHeatmap(SO4, features = top10$gene): The following features were
## omitted as they were not found in the scale.data slot for the SCT assay: H2-Q6,
## Gbp3, Ifi47, Tmem265